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ABSTRACT 

To better understand the nature of the multiphase material found in outflowing galaxies, we study 
the evolution of cold clouds embedded in flows of hot and fast material. Using a suite of adaptive- 
mesh refinement simulations that include radiative cooling, we investigate both cloud mass loss and 
cloud acceleration under the full range of conditions observed in galaxy outflows. The simulations are 
designed to track the cloud center of mass, enabling us to study the cloud evolution at long disruption 
times. For supersonic flows, a Mach cone forms around the cloud, which damps the Kelvin-Helmholtz 
instability but also establishes a streamwise pressure gradient that stretches the cloud apart. If time 
is expressed in units of the cloud crushing time, both the cloud lifetime and the cloud acceleration rate 
are independent of cloud radius, and we find simple scalings for these quantities as a function of the 
Mach number of the external medium. A resolution study suggests that our simulations have sufficient 
resolution to accurately describe the evolution of cold clouds in the absence of thermal conduction 
and magnetic fields, physical processes whose roles will be studied in forthcoming papers. 


1. INTRODUCTION 

Galaxy outflows occur in rapidly star forming galaxies 
of all masses and at all redshifts, and they play a cen¬ 
tral role in the history of galaxy formation {e.g. Heck¬ 
man 1990; Bomans et al 1997; Franx et al 1997; Martin 
1999; Pettini et al 2001; Frye, Broadhurst, & Benitez 
2002; Rupke et al 2005; Veilleux et al 2005; Weiner et al 
2009; Martin et al 2013). They are thought to cause the 
strong correlation between mass and metallicity observed 
in low-mass galaxies {e.g. Dekel & Silk 1986; Tremonti et 
al 2004; Erb et al 2006; Kewley & Ellison 2008); they are 
needed to reconcile the number density of observed galax¬ 
ies with the favored cosmological model (e.g. Somerville 
& Primack 1999; Cole et al 2000; Scannapieco et al 2002; 
Benson et al 2003); and they are essential to the en¬ 
richment of the intergalactic medium (Tytler et al 1995; 
Songaila & Cowie 1996; Rauch et al 1997; Simcoe et al 
2002; Pichon et al 2003; Schaye et al 2003; Eerrara et 
al 2005; Adelberger et al 2005; 2006; Steidel et al 2010; 
Martin et al 2010). Yet, despite the central importance 
of galaxy outflows, the processes that control their evo¬ 
lution are extremely difficult to constrain both theoreti¬ 
cally and observationally. 

Erom a theoretical point of view, the most uncertain is¬ 
sue is the coupling of stars to the surrounding interstellar 
medium (ISM). This is because the highly efficient cool¬ 
ing within the ISM makes it impossible to model super¬ 
novae by adding thermal energy to the medium, while, at 
the same time, the range of physical scales involved does 
not allow for the direct modeling of supernovae within a 
galaxy-scale simulation. As a result, studies have been 
forced to adopt a number of unsatisfactory approxima¬ 
tions, including: temporarily lowering the densities of 
heated particles or delaying their cooling {e.g. Gerritsen 
& Icke 1997; Thacker & Couchman 2000; Stinson et al 

^ School of Earth and Space Exploration, Arizona State Uni¬ 
versity, P.O. Box 871404, Tempe, AZ, 85287-1404 

^ Universitat Hamburg, Hamburger Sternwarte, Gojen- 
bergsweg 112, 21029, Hamburg, Germany 


2006), imposing a minimum temperature floor (Suchkov 
et al 1994; Tenorio-Tagle & Muoz-Tuhon 1998; Strick¬ 
land & Stevens 2000; Eujita et al 2004), using an empir¬ 
ical heating function to mitigate cooling (Mac Low et al 
1989; Mac Low & Eerrara 1999), implementing exagger¬ 
ated momentum kicks (Navarro & White 1993; Mihos & 
Hernquist 1994; Scannapieco et al 2001), and temporar¬ 
ily decoupling particles from their neighbors (Springel & 
Hernquist 2003; Scannapieco et al 2006; Dalla Vecchia & 
Schaye 2008). In fact, even the most detailed cosmologi¬ 
cal ‘zoom-in’ simulations {e.g. Gnedin et al 2009; Agertz 
et al 2009; Geverino et al 2010; Governato et al 2010; 
Shen et al 2012) are faced with the problem that exces¬ 
sive cooling is exacerbated by the fact that supernovae of¬ 
ten go off within giant molecular clouds that must be pre¬ 
conditioned by ionization fronts {e.g. Matzner 2002) and 
radiation pressure {e.g. Murray et al 2010) to be mod¬ 
eled accurately. Although recent efforts have attempted 
to couple stars with the ISM in more realistic detail {e.g. 
Hopkins et al 2011, 2012a,b), even these rely strongly on 
tuning against observations to achieve good results. 

Erom the observational point of view, the most impor¬ 
tant issue is measuring and interpreting the evolution 
of the wide range of multiphase material found in galaxy 
outflows. This ranges from ^ 10^ —lO^K plasma observed 
in X-rays {e.g. Martin et al 1999; Strickland & Heckman 
2007, 2009), to ^ lO^K material observed at optical and 
near UV wavelengths {e.g. Pettini 2001; Tremonti et al 
2007; Martin 2012; Soto & Martin 2012), to 10 — 10^ K 
molecular gas observed at radio wavelengths {e.g. Walter 
2002; Sturm 2011; Bolatto et al 2013). Eurthermore, the 
easiest phase to interpret, the X-ray emitting medium, is 
the most difficult to observe. In fact, the 10^ K medium 
is so hot that it is only detectable in deep Ghandra and 
XMM imaging in very nearby galaxies {e.g. Strickland 
& Heckman 2007; Wang et al 2014) and best measured 
in M82, where it appears to be well fit by a simple an¬ 
alytic model (Ghevalier & Glegg 1985; Heckman et al 
1990). Below 10^ K, the X-ray emitting medium is de- 
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tectable in a larger number of galaxies, and it appears 
to be best understood in terms of ambient material that 
is shock heated by the wind fluid as it flows out to large 
distances {e.g. Suchkov et al 1994; Strickland & Stevens 
2000). Given the high temperatures of the X-ray emit¬ 
ting plasma, it will always escape the gravitational poten¬ 
tial of its host, and its temperature, surface brightness, 
and other properties are remarkably regular over a wide 
variety of galaxies (Grimes et al 2005). 

The colder phases, on the other hand, can be easily ob¬ 
served from the ground and studied at many redshifts, 
but they are poorly understood theoretically. Unlike the 
X-ray emitting medium, low-ionization state material is 
observed to have complex velocity profiles {e.g. Westmo- 
quette et al 2009; 2012) that are strongly correlated with 
the overall host luminosity {e.g. Martin 2012) circular 
velocity {e.g. Martin 2005), star formation rate per unit 
area {e.g. Ghen et al 2010) and star formation rate per 
unit stellar mass {e.g. Heckman 2014). However, because 
it is often visible only through absorption lines and reso¬ 
nant Lyman-alpha emission (Pettini et al 2001), the total 
mass in this phase is poorly constrained. It could either 
be the primary avenue for baryon and metal ejection, or 
make up only a small fraction of the ejected material. 
Equally poorly understood is the final fate of this ma¬ 
terial, as its position along the line of sight is unknown 
and it is often moving at velocities that are similar to 
the escape velocity of the host. In fact, even the pres¬ 
ence of this medium is surprising, as simple theoretical 
estimates predict that it should be disrupted by interac¬ 
tions with the hot wind fluid well before it is accelerated 
to significant velocities {e.g. Scannapieco 2013). 

Furthermore, this uncertainty is somewhat surprising 
given that simulations of cold clouds interacting with a 
hot, high-velocity medium have been carried out by many 
authors, with simulations going back almost 50 years. 
Groups have carried out both two and three-dimensional 
simulations: neglecting radiative cooling (Nittmann et 
al 1982; Klein et al 1994), including radiative cooling 
(Woodward 1976; Mellema et al 2002; Fragile et al 2004; 
Melioli et al 2005; Gooper et al 2009; Marinacci et al 
2010, 2011), including both radiative cooling and ther¬ 
mal conduction (Marcolini et al 2005; Orlando et al 2005, 
2006, 2008; Recchi & Hensler 2007), including the im¬ 
pact of magnetic fields (Mac Low et al 1994; Gregori et 
al 1999, 2000; Fragile et al 2005; Orlando et al 2008; 
Shin et al 2008), and including nonequillibrium chem¬ 
istry effects (Kwak et al 2011; Henley et al 2012). Yet, 
despite the usefulness of these studies in elucidating the 
physics of cold clouds in a hot medium, none of them has 
spanned the range of parameters and timescales neces¬ 
sary to study galaxy outflows. 

The key issue in this case is the ability for the hot wind 
to accelerate cold material, before it disrupts it. In the 
definitive study of the non-radiating, hydrodynamic case, 
Klein et al (1994) showed that if hot material moves past 
a cloud at a velocity 'i;hot that is much greater than the 
cloud’s internal sound speed, it will shred the cloud on a 
‘cloud crushing’ timescale 

ice ^ (1) 

'^hot 

where Rcioud is the cloud radius and 'i^hot/Xo^^ i^ 


locity that the resulting shock moves through the cloud, 
with xo the initial density ratio between the cloud and 
the surrounding medium. Subsequent studies showed 
that magnetic fields and radiative cooling can delay this 
disruption by a few cloud crushing times but not prevent 
it {e.g. Mac Low et al 1994; Orlando et al 2005, 2008). 
On the other hand, to accelerate the cloud to the hot 
wind velocity, the impinging mass must be comparable 
to the cloud mass, which as discussed in more detail be- 
low, takes a time taccei « (47r/3)i?Xud/k-^cioud^hot] « 
Xo^^tec- Because the « 10^ K clouds observed in galaxy 
outflows are in rough pressure equilibrium with the sur¬ 
rounding medium (Strickland & Heckman 2007), xo = 
Tiiot/IO^K = 300 — 10^ such that Gccei ^ 20tcc {e.g. Scan¬ 
napieco 2013). Thus it would appear that cold clouds 
never survive to reach velocities comparable to that of 
the hot wind. 

There is, however, one possible caveat. If cloud cool¬ 
ing is efficient and r’hot exceeds the sound speed in the 
hot medium, Cg^hot, a bow shock develops in front of the 
cloud, which protects it from ablation, both by reduc¬ 
ing heating and squeezing it to even higher densities. 
When combined with radiative cooling, these mitigating 
effects allow the cloud to remain intact for many cloud 
crushing times. In fact, the acceleration of clouds when 
^hot = '^hot/cs,hot is large has been simulated by Gooper 
et al (2009) and Kwak et al (2011), but in both studies 
the cloud was accelerated out of the simulation volume 
before a final answer was obtained. 

In this series of papers, we overcome this limitation by 
carrying out adaptive mesh refinement (AMR) simula¬ 
tions in which we automatically shift the frame of the cal¬ 
culation to the center of mass frame of the cloud once per 
cloud crushing time, allowing us to run until the cloud is 
fully disrupted, regardless of the exterior Mach number. 
In this way, we are able to carry out a suite of simu¬ 
lations that spans the parameter space relevant for hot 
material interacting with the cold clouds that provide us 
with the bulk of the observational constraints on galaxy 
outflows. In this first paper, we study both cloud mass 
loss and cloud acceleration accounting for radiative cool¬ 
ing. Future papers will study: (i) the impact of electron 
thermal conduction {e.g. Gowie & McKee 1977), which 
can provide an additional source of cloud heating, lead¬ 
ing to more rapid disruption, and (ii) the impact of mag¬ 
netic fields, which can extend the lifetime of the cloud, 
both by suppressing thermal conduction perpendicular 
to the field lines {e.g. Braginskii 1965; Bogdanovic et al 
2009; Parrish et al 2009) and stabilizing the KH insta¬ 
bility {e.g. Ghandrasekhar 1961; Walker 1981; Jun et al 
1995; Briiggen & Hillebrandt 2001; Dursi & Pfrommer 
2008). 

The structure of this paper is as follows. In Section 2, 
we describe the physics of cold-cloud hot-wind interac¬ 
tions in more detail, outlining the key parameter space. 
In Section 3, we describe our numerical methods includ¬ 
ing changes of frame and selective de-refinement criteria 
that allow us to obtain results out to many cloud crush¬ 
ing times. In Section 4, we present our simulation re¬ 
sults, give fitting formulae for the mass and velocity of 
the clouds as a function of time and exterior conditions, 
and we discuss the underlying physics that leads to these 
scalings. Gonclusions are given in Section 5. 
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2. PHYSICS OF COLD CLOUD DRIVINC BY HOT GALAXY 
OUTFLOWS 

2 . 1 . Hot Wind Properties 

We are interested in clouds impacted by a hot wind 
expanding rapidly from a starbursting galaxy. For most 
galaxies, this hot medium moves many scale heights per 
Myr, while a typical starburst episode lasts for many 
Myrs {e.g. Greggio et al 1998; Forster Schreiber et al 
2003). Thus, the distribution is expected to be well 
approximated by an equilibrium configuration, and this 
seems to be born out observationally in most objects in 
which reliable X-ray analyses can be made {e.g. Heckman 
1990; Heckman et al 1995; Ott et al 2005; Strickland & 
Heckman 2007; Yukita et al 2012). In this case, assuming 
the expansion is spherically symmetric (or equivalently, 
assuming a conical expansion with a fixed opening angle) 
the equations of mass, momentum, and energy conserva¬ 
tion are 


J.2 


(2) 

dr dr 

5 

(3) 

' ^2 7 Lt 1 7 -PhotV 

V 2 7-1 Phot J \ 

= 

(4) 


where phot, ^hot, and 7 are the density, radial ve¬ 
locity, and ratio of specific heats of the hot medium, and 
the mass and energy input rate are 


f gm,o if r < i?* 
\0 ifr >0 


and Qe = hot, 0/(7 - 1)5 respectively, where 

is the driving radius of the flow and C5^hot,o = [(7 ~ 
l)ge,o/^m,o]^^^ is the sound speed of the hot medium at 
r = 0. For reference, in the case of M82, ^ 300 pc 

(Strickland & Heckman 2009). 

The solution to these equations is the Chevalier and 
Clegg (1985) model, which can be used to gain a good 
understanding of the hot wind conditions as a function of 
radius. For ay = 5/3 gas, the sound speed when r < R^ 
is approximately constant at 


Cs,hot (f) 


^s,hot,0 — 0.82 



(6) 


where E is the total energy input from supernovae per 
unit time, e is the fraction of this energy that is deposited 
into the hot medium, rh is the total mass ejected by su¬ 
pernovae per unit time, and (3 is the mass input per unit 
time into the hot material, which may exceed one due to 
entrainement. The outward velocity of the medium, on 
the other hand, increases approximately linearly in the 
driving region as v\iot{r) = Cs^hot,oF/4i?*, such that the 
Mach number increases as Mhot(F) ~ r/AR^^ where here 
and below we use M to distinguish Mach number from 
mass, which is always denoted as m. 

At the edge of the driving region, the solution 
reaches the sonic point, with i;hot(^*) = Cs,hot{Ri.) = 
{eE. Finally, outside the driving region, the 
radial velocity quickly approaches a constant value of 


'^hot(F) = (2e^//dm)^/^ and the Mach number goes 
as Mhot(F) = 2^/^(r/i?*)^/^, such that Cs^\,ot{r) = 
2-5/3( 2 ^^/^p^) 1 / 2 (^/j^^)- 2 / 3 ^ Note that at all radii, the 
Mach number is purely a function oi r/R^ and indepen¬ 
dent of the rate of mass and energy input into the hot 
wind. 

Thus sampling a series of Mach numbers with a suite of 
simulations corresponds to sampling the properties of the 
hot wind as a function of radius. These distances have 
been tabulated in column 3 of Table 1 for a range of 
Mach numbers. Note that the velocity and sound speed 
change rapidly near the sonic point, such that Mhot goes 
from 0.5 to 3.5 between 0.9 and 1.1 r/R^. The Chevalier 
and Clegg (1985) model also gives an estimate of the 
density of the expanding medium, which is indicated in 
column 4 of Table 1 as n/no, the number density of the 
medium relative its number density at r = 0. Again, for 
reference, the pressure of the ~ lO^K medium in M82 is 
P/k ^ 1 — 3 X 10^ cm“^ K in the driving region, which 
corresponds to no ~ 0.1 — 0.3 cm“^. 

If we take a standard estimate of E/rh of 
lO5iergs/(lOM0), this gives a velocity {2eEII3rfi)^^‘^ = 
(e/yd)^/^ 3160 km s“^ at large radii and a temperature 
at small radii of {e/(3) 1.1 x 10^ K or 9.3 keV. This is 
only slightly higher than observed in M82 (Strickland & 
Heckman 2009), implying only moderate mass loading is 
typical for many starbursts, with e/f3, usually > 0.3 and 
always >0.1. This gives us a reasonably small parame¬ 
ter space of conditions of the exterior flow that are most 
important for understanding cold-cloud acceleration. 

Table 1 also shows 'r’hot and Thot for a number of choices 
of e/yd, for the cases described in more detail below. 
These values were chosen to span a wide range of efficien¬ 
cies, focus on the poorly-studied Mhot > 1 regime, and 
provide multiple runs with the same temperature and/or 
velocities, such that we can study the impact of changing 
Mhot while maintaining Vhot or Cs^hot,o constant. Thus we 
only consider a single choice of ej(3 for Mhot = 0.5 case, 
but we consider three choices for Mhot = 1- oae that 
gives a similar 'i;hot as the Mhot = 0.5 case, one that 
gives the same Cg^hot as the Mhot = 0.5 case, and one 
for which e/yd ~ 1. Similarly, for Mhot ~ 3.5 we adopt 
three choices of e/yd that allow for comparisons with the 
Mhot = 1 cases with the same Cs^hot,o and similar Vhot 
values. Note that in one case, this forces us to take e/yd 
slightly greater than 1, which is not unreasonable given 
the uncertainties in the mass and energy input from su¬ 
pernovae. Finally, in the highest three mach number 
cases we chose e/yd to exactly match the Cs7ot,o and e’hot 
taken in the Mhot ~ 3 cases. Together these choices allow 
us to define a small number of runs that can nevertheless 
be used to understand the evolution of cold clouds in a 
wide range of outflow conditions. 

2.2. Cold Cloud Disruption 

The clouds of interest will have temperatures ~ 10^ K, 
due to the balance between photoheating by the ioniz¬ 
ing background and strong radiative cooling between 10^ 
and 10^'^ K. For a medium at lO^K, the Jeans length, 
Aj ~ Cs,cioud(Gp)“^/^ ^ 2kpc(n/cm^)“^/^, meaning 
that the clouds will be pressure-confined rather than 
gravitationally bound. For such clouds, the initial ra- 



4 


TABLE 1 

Simulation Parameters 


Name 

Mhot 

r/R^ 

nino 

A/? 

'^hot 

km 

Aiot 

10® K 

Xhot 

keV 

XO 

-Mcioud 

tec 

Myr/lOOpc 

-^cool 

cm“^ 

^cool 

Mopc-2 

M0.5v430 

0.5 

0.9 

0.8 

0.2 

430 

30 

2.7 

3000 

28 

12.5 

IcF^ 

5 X 10-1 

Mlv480 

1.0 

1.0 

0.4 

0.1 

480 

10 

0.86 

1000 

32 

6.4 

1017 

5 X 10-4 

Mlv860 

1.0 

1.0 

0.4 

0.3 

860 

30 

2.7 

3000 

57 

6.2 

1017 

5 X 10-4 

Mlvl500 

1.0 

1.0 

0.4 

0.9 

1500 

100 

8.6 

10000 

100 

6.5 

1017 

5 X 10-4 

M3.8vl000 

3.8 

1.1 

0.2 

0.1 

1000 

3 

0.27 

300 

66 

1.7 

1017.5 

1.5 X 10-® 

M3.5vl700 

3.5 

1.1 

0.2 

0.4 

1700 

10 

0.86 

1000 

110 

1.8 

1017.5 

1.5 X 10-® 

M3.6v3000 

3.6 

1.1 

0.2 

1.1 

3000 

30 

2.7 

3000 

200 

1.8 

1017.6 

1.5 X 10-3 

M6.5vl700 

6.5 

1.9 

0.05 

0.3 

1700 

3 

0.27 

300 

110 

1.0 

lOi® 

5 X 10-3 

M6.2v3000 

6.2 

1.9 

0.05 

1.0 

3000 

10 

0.86 

1000 

200 

1.0 

lOi® 

5 X 10-3 

M11.4v3000 

11.4 

2.6 

0.03 

0.9 

3000 

3 

0.27 

300 

200 

0.56 

lOi® 

5 X 10-2 


tio of the cloud density to the density of the exterior 
medium, xo, will also be equal to the ratio of the exte¬ 
rior temperature to the cloud temperature. This allows 
us to associate a single density contrast with each exte¬ 
rior temperature, as shown in column 9 of Table 1. The 
fixed temperature of the clouds also allows us to com¬ 
pute a single cloud Mach number Mdoud = '^hot/cs,cloud? 
shown in column 10 of this table. 

As the sound crossing time for these pressure con¬ 
fined clouds is only ^ 0.05 Myr i^cioud/parsec, we can 
reasonably approximate them as spherical before they 
encounter the shock. This assumption also makes the 
simulations of their evolution easier to interpret, as well 
as allows us to make closer contact with the rich exist¬ 
ing literature on the topic. There are a few important 
timescales that determine the evolution of cold clouds 
in this case. The most important of these is the cloud 

crushing time mentioned above, tec = ^cioud/('^hotXy^)- 
In the Appendix, we relate Uhot? the velocity of the wind 
impacting the cloud to Ut, the velocity of the initial shock 
transmitted through the cloud, and show that in the limit 

of high Mach numbers Vt ~ '^hot/Xo^^- Thus the cloud 
crushing time gives a rough estimate of when the shock 
from the exterior medium moves through the cloud, heat¬ 
ing it and disrupting it if it is not able to radiate the en¬ 
ergy away efficiently. In column 11 of Table 11 we show 
the value of this time in units of Myrs for clouds of size 
-^cloud ~ lOOpC. 

A second important timescale is the cooling time be¬ 
hind the transmitted shock, which can be estimated as 

^ _ 3/2ncioud^Tt 

tcool x /rri ^ ’ V'/ 

t j'^e,cloud'll,cloud 

where k is the Boltzmann constant, Tt is the post-shock 
temperature, A(Tt) is the equilibrium cooling function 
evaluated at Xf, and ricioud? ^e,cloud? and n^^cioud? ^^re 
the total, electron, and ion number densities within the 
cloud, respectively. The ratio of this time to the cloud 
crushing time is tcooi/icc = -Ncooi/(«i,cioudi'c), where the 
column density, A^cooi = 3A:Ttuncioud[2Axi/2ne,cioud]“^ is 
purely a function of the velocity of the transmitted shock. 

Because the clouds of interest are at a fixed temper¬ 
ature and the Mach number of the transmitted shock 
within the cloud can be related to the Mach number 
in the exterior medium, this means that A^cooi can be 
well approximated by a function of Mhot- These num¬ 
bers have all been added to Table 1, as calculated using 


the equilibrium cooling curves in Wiersma et al (2009), 
assuming solar metallicity and an average atomic mass 
of 0.6. These values indicate that for the conditions of 
interest, the cooling time will be much shorter that the 
cloud crushing time for all but the smallest, sub-parsec 
size clouds, and thus radiative cooling will be very effi¬ 
cient throughout the evolution of most observed clouds. 

A third important timescale is the characteristic 
timescale for disruption by the shear-driven Kelvin- 
Helmholtz (KH) instability. In a linear stability anal¬ 
ysis, the width of the KH layer grows at a rate Aukh oc 
where Au is the velocity difference between the 
cloud and the hot medium (Chandresekshar 1961). Sim¬ 
ilarly, experimental measurements of subsonic shearing 
flows show that the nonlinear growth across the bound¬ 
ary with a high density ratio is asymmetric, such that 
the “entrainment ratio” of the width on the high density 
side of the layer to the width on the low density side of 
the layer Ey = x~^^^ (Brown 1974; Konrad 1976). 

This can be understood in terms of entrainment into a 
spatially-growing shear layer made up of large-scale vor¬ 
tical structures convecting at a velocity Vc (Coles 1981; 
Dimotakis 1986). In the subsonic case, Vc is set by the 
condition that in the frame of the vortices 

Pl+pl{Vi-Vcf K p2{V2-Vc)‘^ + P 2 , ( 8 ) 

where Pi and P 2 are the pressures of the material on 
either side of the layer, pi and p 2 are their densities and 
vi and V 2 are their velocities. This means that if Pi = P 2 ? 
then \vi—Vc\ = X^^^l'^ 2 ~'^c|? such that the heavy material 
drags the vortices along with it. The expansion of the 
shear layer on either side of the density contrast in this 
case appears to be roughly proportional to \v—Vc\ (Brown 
& Roshko 1974; Papamoschou & Roshko 1988; Slessor et 
al 2000), leading to a strongly asymmetric growth of the 
boundary with the width on the dense side as a function 
of time given by 

S ^ O.ltAuy”^/^. (9) 

Thus, even if cooling is efficient, the subsonic KH insta¬ 
bility will grow to the scale of the cloud within a few 
cloud crushing times. 

In the supersonic case, however, the growth of the 
boundary layer is much slower (Chinzei et al 1986; 
Papamoschou & Roshko 1988; Samini & Elliott 1990; 
Goebel & Duttonl991; Hall et al 1993; Barre et al 1994; 
Clemens & Mungal 1992; Naughton et al 1997; Slessor 
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et al 2000), but the level of this suppression, and the 
structure of the mixing layer are much more poorly un¬ 
derstood. Papamoschou & Roshko (1988) measured this 
suppression as a function of convective Mach numbers: 


Me,i = 


Vl - Vc 
Cs,l 


and Mc 2 = — —(10) 

Cs,2 


where again Vc is the convection velocity of the turbu¬ 
lent structures and and Cs ^2 are the sound speeds of 
the two media. Although in this case the best choice 
for Vc is unclear. If one adopts the isentropic pressure- 
recovery model of Papamoscho & Roshko (19881, then 
if 7 i = 72 , Vc is the same as given by eq. ^ and 
Mcp = Me ,2 = Au/(cs,i -h Cs, 2 )- As a function of this 
quantity, the width of the layer is reduced by a factor of 
6{Mc)/S{0) which approaches 1/4, as Me approaches the 
largest experimentally measured values ^1.5. 

On the other hand, several authors have argued for 
different choices for Uc, such as the phase speed of the 
linearly most unstable mode (Ragab & Wu 1989), the 
phase speed of the linearly most unstable mode at a point 
corresponding to neutral stability (Sandham & Reynolds 
1989), and the speed of turbulent structures in the pres¬ 
ence of asymmetric shocks (Dimotakis 1991). In fact, 
experimentally it appears that supersonic shear layers 
favor a configuration in which Mc,i and Me ,2 are very 
different, such only one side of the layer is supersonic 
(Papamouschou 1991). With this asymmetry in mind, 
Slessor et al (2000) proposed a suppression of the width 
of the KH layer by a factor 

S{M,)/5{0) « [1 + 4(7 - (11) 

where Mkh = Au/cs,min is the ratio of Av with Cs,min 5 
the minimum of the two sound speeds of the fluids 
that are being sheared (Slessor et al 2000). Extrapo¬ 
lating this fit formula to the cases of interest reduces 
the growth of the shear layer by an even large factor 
than the Papamoschou & Roshko (1988) scaling, with 
^(Mkh)/^( 0) ~ Mg,cloud at very large Mach numbers. 


2.3. Velocity Evolution 

The radial velocity of the cloud, Udoud, will depend on 
the momentum imparted by the hot wind, divided by the 
total cloud mass. For a spherical cloud this gives 


'^cloud(^) — / 

Jo 


+f\2 


/7rRcioud(^') Au(p) 

fi?cioud(i')"x(iO 


3uhot 


ptltc 

Jo 


di' R\{i')Av{i 


'f\2 


( 12 ) 


where i is the time in units of the initial cloud crushing 
time, xo is again the initial density contrast, and Av and 
R± are the relative velocity between the cloud of the ex¬ 
terior medium and the radius of the cloud perpendicular 
to the flow in units of their initial values, respectively. 
Notice that because this ratio compares the cross section 
the cloud presents to the incoming flow with the total 
mass of the cloud, the final expression does not depend 
on the size of the cloud in the direction of the incoming 
material. Notice also that if we consider the evolution of 
the cloud in units of the cloud crushing time, the physi¬ 


cal cloud radius does not appear in this expression, only 
~ R-L / Rc\ond • 

Eq. ( [T^ also illustrates the fact that it is more dif¬ 
ficult to accelerate clouds with larger density contrasts, 
although this effect is somewhat mitigated because the 
cloud is initially flattened by the collision {e.g. Klein et 
al 1984). In the limit in which Ah, and R± are fixed, 
the distance traveled by the cloud is proportional to the 
square of its lifetime in units of the cloud crushing time. 
This makes the results especially sensitive to the late 
time evolution of the cloud, and it also means that the 
suppression of the KH instability can be very important 
in determining the distance the cloud travels before it is 
disrupted. 


3. METHODS 
3.1. Setup 

To study cold-cloud hot-flow interactions over the 
range of conditions encountered in galaxy outflows, we 
carried out a suite of simulations using FLASH (version 
4.2), a multidimensional hydrodynamics code (Fryxell et 
al 2000) that solves the Riemann problem on a Cartesian 
grid. Each of these runs adopted one of the sets of repre¬ 
sentative conditions described in Table 1. All simulations 
were three-dimensional, as imposing a 2D cylindrical ge¬ 
ometry reduces the degrees of freedom over which shear- 
driven instabilities can develop, which is likely to have a 
large impact on the evolution of the cloud {e.g. Pan et al 
2012 ). 

In all cases, we used the default directionally-split 
(Strang 1968) Piecewise-Parabolic Method hydrody¬ 
namic solver (PPM; Colella & Woodward 1984; Colella 
& Glaz 1985; Fryxell, Muller, & Arnett 1989). We chose 
this approach as we found that it was much better at pre¬ 
serving spherical symmetry in test simulations in which 
a cloud was compressed by a stationary high-pressure 
medium. We also made use of the shock detect flag 
which lowered the prefactor in the Courant-Friedrichs- 
Lewy timestep condition from its default value of 0.4 to 
0.25 in the presence of strong shocks. 

Because of the scalings discussed in §2, as long as the 
column density is sufficiently large that cooling is effi¬ 
cient, the size of the cloud scales out of the problem 
if we express our results in units of the cloud crush¬ 
ing time. Thus without a loss of generality, we chose a 
fixed cloud radius of 100 parsecs for all our simulations, 
an initial temperature of lO^K, and a mean density of 
p = 10“^^ g cm“^, such that pRcioud = 1-4M0 pc“^, 
and ni,cloud^cioud = 3.1 x 10^^ cm“^. In all simulations, 
the computational domain covered a physical volume of 
-800 X 800 parsecs in the x and y and directions and 
—400 to 800 parsecs in 2 : direction, where the cloud was 
initially centered at (0,0,0) and ^ is the direction of the 
hot outflowing material. Outside of the cloud, the initial 
velocity and the sound speed of the material were taken 
to be Uhot: and Cg^hot as given by Table 1, and the density 
was set by pressure equilibrium with the cold cloud. 

At the lower 2 : boundary we continuously added mate¬ 
rial to the grid with the same values of Uhot^ and 

density as in our initial conditions. In the x and y di¬ 
rections, as well as at the -hz boundary, we adopted the 
FLASH “diode” boundary condition, which assumes a 
zero normal gradient for all flow variables except pres- 
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sure and does not allow material to flow back onto the 
grid. The large x and y size of the simulation was chosen 
such that this zero gradient condition did not effect the 
shape of the shock, which can become unnaturally pla¬ 
nar in smaller simulation domains. Likewise a significant 
standoff distance was left in front of the cloud, such that 
the front of the bow shock remained within the simula¬ 
tion volume, and a significant distance was left behind 
the cloud, to capture the evolution of the disrupted ma¬ 
terial. 

3.2. Cooling, Frame Changing, and 
Refinement/D e-refinement Criteria 

Cooling was computed in the optically-thin limit, as¬ 
suming local thermodynamic equilibrium as 

= (13) 

where ^cooi is the radiated energy per unit mass, p is the 
density in the cell, rUp is the proton mass, Y = 0.24 is the 
helium mass fraction, /i = 0.6 the mean atomic mass, and 
A(T, Z) is the cooling rate as a function of temperature 
and metallicity. Here we made use of the tables compiled 
by Wiersma et al (2009) from the CLOUDY code (Fer- 
land et al 1998), making the simplifying approximations 
that the metallicity of the material is always solar and 
that the abundance ratios of the metals always occurs 
in solar proportions. As in Gray & Scannapieco (2010), 
subcycling was implemented within the cooling routine 
itself, such that T and A(T, Z) were recalculated every 
time Ecoo\/E > 0.1. This is equivalent to an integral 
formalism that assumes a constant density over each hy¬ 
drodynamic time step {e.g. Thomas & Couchman 1992; 
Scannapieco, Thacker, & Davis 2001). We did not in¬ 
clude heating by a photoionizing background in our cal¬ 
culations, but the approximate impact of including this 
would be to slightly raise the minimum temperature be¬ 
low which the cloud can cool effectively, while having a 
minor effect at higher temperatures. Likewise changing 
the metallicity would lengthen the cooling times some¬ 
what, particularly above 10^‘^K, but not change the over¬ 
all conclusion of the efficiency of cooling for the majority 
of the clouds of interest. 

In order to be able to run our simulations until the 
clouds were disrupted, while still keeping them on the 
grid, we developed an automated frame changing rou¬ 
tine. By labeling the cloud material with a scalar, we 
were able to track its center of mass position and veloc¬ 
ity, Xcioud and Vcioud, as well as its radial extent in the 
X, y, and z directions, calculated as the mass weighted 
average values of abs{x — Xdoud), abs(^ — Xcioud), and 
ahs{z — ^^cioud)- Every cloud crushing time, starting at 
t = 3tcc, we checked if the cloud was moving in the pos¬ 
itive z direction in the frame of the simulation and if 
^cioud ^ ^cioud/2, such that the cloud was not too near 
the front of the simulation volume. If these criteria were 
satisfied, we shifted every z velocity in the simulation by 
a constant factor of 1.2 times the instantaneous cloud ve¬ 
locity, appropriately adjusting the kinetic energy. This 
moved the simulation from a frame in which the cloud 
was rapidly drifting towards the -\-z boundary, to a new 
frame in which the cloud was slowly drifting towards the 
—z boundary. As more momentum accumulated on the 


front of the cloud, it would once again start to move 
toward the boundary, until another cloud crushing 
time went by and the frame was again shifted to make 
the cloud slowly drift upstream. In this way, we were 
able to keep the cloud center of mass near z = 0 at all 
times, and still reconstruct its evolution in the frame of 
the wind by keeping track of the overall frame shift. 

To determine when zones are refined and de-refined, 
FLASH uses the second derivatives of “refinement vari¬ 
ables,” normalized by their average gradient over a cell. 
In the default FLASH configuration, if this number is 
greater than 0.8, the cell is marked for refinement, and if 
all the cells in a block lie below 0.2, they are marked for 
de-refinement. We used both density and temperature as 
refinement variables in this way, but also adopted a set 
of additional refinement and de-refinement criteria, cho¬ 
sen to minimize the computational cost of the simulation 
while at the same time maintaining the most accurate 
results possible in the spatial regions that are the most 
important to the evolution of the cold cloud. 

While the large simulation volume perpendicular to 
the z axis was needed to maintain the proper angle for 
the Mach cone that forms around the cloud, it is not 
important for us to achieve high resolution for this pur¬ 
pose. Thus we forced our simulation to automatically 
mark cells for de-refinement if either of two criteria were 
satisfied: (i) if the distance from a cell to the z axis was 
greater than three times the original cloud radius or nine 
times the instantaneous x extent of the cloud or if (ii) 
the distance from a cell to the z axis was greater than 
either the original cloud radius or three times the instan¬ 
taneous X extent of the cloud and also both abs(z) and 
ahs{z — ^^cioud) were greater than 3 i^cioud- 

Likewise, in order to obtain the most reliable results 
possible, it is important for us to maintain high reso¬ 
lution consistently in the regions within the cloud and 
the shear layer immediately around it. To achieve this 
in the initial stages of the interaction, when t < 2tcc, 
we automatically marked cells for maximum refinement 
within the cylindrical region with abs(z) < 1.5i?cioud and 
the distance from the z axis {x‘^ + < l-5i?cioud- 

Once the transmitted shock made its way all the way 
through the cloud at t > 2tcc, we then only forced 
the simulation to maintain maximum refinement when 
abs(2: ^cioud) ^ -^cioud and (x T ^ ) ^ l*6-Rcioud5 

although, in practice the complex density and pressure 
structure that developed around the cloud by this time 
led to highly-refined regions that occupied a much larger 
volume, approaching the full volume in which derefine¬ 
ment was not forced by our criteria above. 

4. RESULTS 

We carried out twelve simulations in all: ten runs with 
the same maximum resolution, which were used to span 
the parameter space in Table 1, and a two additional runs 
with different maximum resolutions, which were used to 
examine resolution effects. Each of our parameter study 
runs was carried out on a base 64 x 64 x 48 grid, with four 
additional levels of refinement. The base resolution for 
these runs was 25 parsecs in each direction and the max¬ 
imum resolution was 1.5625 parsecs, or i?cioud/64. This 
value was chosen to be significantly smaller than the size 
advocated by MacLow & Zahnle (1994), who found that 
their 2D cylindrical simulations of the breakup of comet 
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Fig. 1.— Mass evolution of the cloud as a function of mach number and velocity. In each panel the solid blue line shows the 

fraction of the mass at or above 1/3 the original cloud density, the dotted red line shows Fi, the fraction mass at or above the original 
cloud density, and the purple dashed line shows the fraction of the mass at or above 1/10 of the original density. In all panels the 

time is given in units of the cloud crushing time, and the length of each panel is proportional to the time in Fc units. 


Shoemaker-Levy 9 in Jupiter’s atmosphere required at 
least 25 zones across the projectile to achieve convergent 
results. At this resolution, our simulations took 10-45k 
CPU hours to run, depending on the choice of parame¬ 
ters. 

4.1. Mass Evolution 

In Figured we show the fraction of the mass retained 
by the cloud^as a function of time in units of the cloud 
crushing time. In this figure, the retained mass frac¬ 


tion is defined using three different measures: Ui(t), the 
fraction of the total mass at or above the original cloud 
density, Fi/3(t), the fraction of total mass at or above 
1/3 the original cloud density, and Ui/io(t), the fraction 
of the total mass at or above 1/10 of the original cloud 
density. Although Fi is the simplest measure of cloud 
mass, it displays a number of misleading features. Our 
simulations start with simple initial conditions, with the 
exterior medium moving at all points outside the cloud. 
Thus the back of the cloud expands for roughly 0.5tcc, 






































while the flow behind the cloud rearranges itself. This 
causes all runs to display an initial drop in Fi as a sub¬ 
stantial part of the cloud drops slightly in density. Con¬ 
versely, from about tec to 2 tcc, the shock at the front of 
the cloud increases the cloud density. In the high Mach 
number runs, the density increase behind this (radiative) 
shock is large, bringing almost the entire cloud back to a 
density equal to or exceeding the original density, which 
brings Fi back to ~ 1 by t = 2tcc- At lower Mach num¬ 
bers, the density increase is somewhat smaller, such that 
the increase in Fi is more modest. 

These transient features are avoided by working with 
Fi /3 or Ti/iQ. In these cases, the initial downstream ex¬ 
pansion and shock compression do not cause a significant 
change in the material interpreted as belonging to the 
cloud. Instead, the cloud evolution is largely monotonic, 
remaining approximately constant for the initial stages of 
the interaction, and then dropping at a roughly constant 
rate at times ^ 2tcc- This evolution is extremely simi¬ 
lar between F 1/3 and Ti/ 10 , indicating that the precise 
choice of density threshold does not affect our conclu¬ 
sions. On the other hand, the timescale on which these 
mass fractions evolve varies extensively between the runs, 
such that the cloud in the Mhot = 0-5, v = 430 km/s run 
has lost 3/4 of its mass by t = 5tcc, but at this same 
time, the cloud in the Mhot = llA^ v = 3000 km/s case 
still contains over 90% of its original mass. 

To better illustrate the evolution of the clouds, in Fig¬ 
ures ii and 1^ we plot slices of the central density dis¬ 
tribution at four characteristics times: tgo, the time at 
which Fi /3 = 90%, ^ 75 , the time at which F 1/3 = 75%, 
tso, the time at which F 1/3 = 50%, and ^ 25 , the time at 
which Fi/3 = 25%. In the lowest Mach number runs, il¬ 
lustrated in Figure tgo usually occurs during the early 
stages of the interaction, when t ^ 2tcc- For example, in 
the Mhot = 0.5, 1 ; = 430 km/s run, the cloud loses 10% of 
its mass by t = 1 . 6 tec, mostly due to material stripped 
from the sides during the initial interaction. As time 
goes on, this run evolves primarily through the growth 
of the Kelvin-Helmholtz instability, which works its way 

in from the sides at a rate oc 'T’hot/Xo^^’ disrupting the 
majority of the cloud by 5tcc- 

As expected from the discussion in § 2 . 2 , the growth of 
the KH instability is somewhat slowed in the Mhot = 1 
case, such that ^25 ~ Stec with slightly higher values be¬ 
ing seen in the highest velocity cases. This evolution is 
extended even further in the Mach 3.5 runs, which 
retain over 90% of their mass until well after the ini¬ 
tial stages of the interaction. Thus the slices taken at 
tgo in Figure show the clouds in highly shocked and 
compressed configurations, with numerous features that 
have densities well above the initial density of p = 10 “^^ 
g cm“^. As time evolves in the supersonic runs, these 
clouds become more cometary in appearance, with the 
densest features located near the head of the cloud, which 
is trailed by a long tail of lower density material. This 
tail is then slowly mixed into the surrounding medium 
such that most of the mass is lost by lOtcc- 

These features are exaggerated even further in the 
Mach « 6.5 and 11 cases, shown in Figure Here the 
tgo plots, taken at ~ 5tcc, show regions with densities ex¬ 
ceeding 30 times the original density. Again these clouds 
evolve into cometary distributions with moderate den¬ 


sity contrasts, and these retain a large fraction of their 
mass for even longer than the Mhot ~ 3.5 runs, such that 
^25 ~ 20tcc- Note that these long disruption times occur 
at t/tec values later than have been simulated previously. 
On the other hand, they occur sooner than expected from 
the laboratory measurements, especially as fit by Slessor 
et al ( 2000 ) with eq. ( 11 ). 

In Figure we show tgg, ^75, tso, and ^25 in units of 
cloud crushing times for each of our runs, plotted as a 
function of Mach number. In all cases, we find that the 
scatter between runs with difference xo values but the 
same Mach number is small, meaning that the differences 
in disruption times at constant Mhot is well accounted for 

by the scaling of tec- Similarly, the scaling with Mach 

number is well fit by a simple function, such that in all 
panels we obtain good agre ement wi th our simulations 
results and timescales oc + Mhot- In particular, we 
find that our results match 


t — atec y/l + Mhot 


(14) 


where a= 1.75, 2.5, 4, and 6 , at tgg, t75, tso, and t25, 
respectively. 

To better understand the origin of this scaling, in Fig¬ 
ure we plot slices of the pressure, density, and velocity 
distribution from runs with four different Mach numbers: 
an Mhot = 1 run with Vhot = 840 km/s, and three runs 
with Vhot = 3000 km/s and Mhot = 3.6, 6.2, and 11.4. 
For ease of comparison between runs, all slices in this 
figure are plotted at the same point in their evolution, 
when t/tec = 6 . From the top panels in this figure, we can 
directly compare the pressure distribution established in 
each of these runs. Here we see that the increase in pres¬ 
sure is strongest at the front of the cloud, with values 
approaching the ^ 1 + M^ increase expected for a nor¬ 
mal shock [see eq. (A 2 )]. 

However, downstream from the front of the cloud the 
pressure increase is more modest, scaling instead as ~ 
1 + M. This weaker scaling with Mach number can be 
understood by the fact that the majority of the cloud 
material does not encounter flow material that has passed 
through a normal shock, but rather through an oblique 
shock. In the case in which flow is passing over a wedge 
with opening angle an oblique shock is formed with a 
downstream pressure given by 


^ - 1 + 

Pi ~ \/M2 - 1 


1 + 76 »M 


(15) 


(Ackert 1925), rather than While this is a some- 

what idealized example, the strong alignment between 
the shock front and the direction of the hot material 
suggests that the flow is much better approximated by 
an oblique shock encountering a 6 > I /7 30° wedge 

than a normal shock encountering a ^ ~ 90° wall. This 
large difference in pressure between oblique and normal 
shocks also highlights the importance of adopting a large 
simulation domain, such that the zero normal gradient 
assumed for all flow variable at the boundaries does not 
affect the shape of the shock. 

By comparing the pressure slices to the density slices, 
we can also see that the downstream shock is well- 
separated from the edge of the cloud. The primary role 
of the shock throughout the majority of the cloud is to 
increase its density by 1 + M, as this material must 
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Fig. 2.— Plots of central density slices from simulations with Mach numbers 0.5 and 1, at times at which the fraction of the mass at 
or above 1/3 the original density of the cloud is 90% (tgo, first column), 75% (^ 75 , second column), 50% (tso, third column), and 25% 
(^25, fourth column). First row: Results from simulation M0.5v430, at times tgo = l.Otcc, ^75 = 2.2tcc, ^50 = 3.2tcc, and ^25 = 4.8tcc- 
Second row: Results from simulation Mlv480 at times tgo = l-8tcc, ^75 = 2.6tcc, ^50 = 4.8tcc, and ^25 = T.Otcc- Third row: Results from 
simulation Mlv860 at times tgo = 2.0tcc, tr^ = 3.4tcc, tso = O.Otcc, and t 25 = 7.6tcc- Fourth row: Results from simulation Mlvl500 at 
times tgo = 2.4tcc, tr^ = 3.8tcc, tso = 5.8tcc, and t 25 = 8 . 8 tcc- All densities are in g cm“^ and all lengths are in kpc. These figures, along 
with the slice plots below, are insets of the region 0.6 x 0.8 kpc region around the cloud, which is a subset of the full 1.6 x 1.6 x 1.2 kpc 
simulation volume. 
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Fig. 3.— Plots of central density slices from simulations with Mach numbers 3.6. As in Figure the results are shown at times at 
which the fraction of the mass at or above 1/3 the original density of the cloud is 90% (tgo, first column), 75% (trs, second column), 50% 
(tso, third column), and 25% (^25, fourth column). First row: Results from simulation M3.8vl000 at times tgo = 5.0tcc, ^75 = 6 . 6 tcc, 
^50 = 8-4tcc, and ^25 = 11.6tcc- Second row: Results from simulation M3.5vl700 at times tgo = 4.2tcc, ^75 = 5.8tcc, tso = 7.6tcc, and 
t25 = 10.4tcc- Third row: Results from simulation M3.6v3000 at times tgo = 4.2tcc, tr^ = G.Otcc, t^o = S.Otcc, and t 25 = 12.0tcc- All 
densities are in g cm“^ and all lengths are in kpc. 
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Fig. 4. — Plots of central density slices from simulations with Mach numbers 6.5 and 11, with times as in Figures and 
First row: Results from simulation M6.5vl700 at times tgo = 4.4tcc, ^75 = 6.4tcc, tso = lO.Gtcc, and ^25 = IT.Otcc- Second row: 
Results from simulation M6.2v3000 at times tgo = 4.4tcc, ^75 = 6.4tcc, tso = 9.6tcc, and ^25 = 19.0tcc- Third row: Results from simu¬ 
lation M11.4v3000 at times tgo = 5.6tcc, trs = 8.8tcc, too = 13.0tcc, and t 25 = 22.2tcc- All densities are in g cm“^ and all lengths are in kpc. 


rearrange itself to achieve pressure equilibrium with the 
exterior flow while still maintaining a constant temper¬ 
ature. Thus a large pressure gradient develops between 
the head and tail of the cloud, which leads to an accel¬ 
erating expansion of the cloud in the 2 ; direction by a 
rate 


Avz{t) 


t dP _ ^ Photvlot 

Pcloud Pcloud-^cloud 


such that 


t l^hot 


G (1 + MhotlXo^^ 
(16) 


Rzjt) f ^ y 

^cloud 2 \tcc\/l + ^hot / 


This stretching in the direction of the flow is the source 
of the cloud’s cometary appearance at late times in the 
high Mhot runs. Furthermore, from our empirical fit in 


eq. (14) we find that Rz ~ is almost purely a function of 
the mass lost from the cloud, such that Rz ~ l-Si^cioud, 
3 -Rcioud 5 S-^cioud? and 18-Rcioud ^90 5 ^75 5 ^50 5 and ^ 25 , 
respectively. 

These estimates are roughly consistent with the lengths 
seen in the density slices, with the exception of the dis¬ 
tributions at ^ 25 , in which it is very difficult to assign 
lengths to the clouds, which have become fragmented 
into many small clumps. Thus, even as the KH in¬ 
stability is strongly damped in these supersonic runs, 
the clouds are also subjected to rapid stretching by the 
streamwise pressure gradient. This first distorts the 
clouds, then finally splits them apart as Rz{t)/Rcioud 
grows rapidly on a timescale (x tccVl + Tfhot- 


4.2. Velocity Evolution 

Next we consider the velocity evolution of the cloud, 
comparing our simulation results to the simple seal- 
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Fig. 5.— Time at which the fraction of the mass at or above 1/3 
the original density of the cloud is 90% (tgo), 75% (trs), 50% (tso), 
25% (^ 25 ) in units of the cloud crushing time. The (violet) points 
show the simulation results, and the (blue) lines show the fits tgo ^ 
1.75tccVl + -^h ot; ^75 ^ 2.5tccVl + -^hot 7 ^50 ~ ^tccV^ + ^hot 5 
and t25 ~ 6tcc\/l + ^hot? respectively. 


ings expected from eq. ( 12 ) above. Because the pri¬ 


mary impact of the shocks seen at high Mhot values is 
to compress the cloud, this means that the extent of 
the cloud perpendicular to the direction of the shock, 
R±{t) = R±{t)/Rcioud^ will drop as the cloud is squeezed 
by the pressure increase behind it. Thus, although the 
high Mach number cases are able to avoid disruption for 
many more cloud crushing times, the acceleration of such 
clouds per cloud crushing time will also be smaller. In 
Figure]^ we plot the velocity evolution of each our sim¬ 
ulations, using units of both km/s and Uhot- Here we 
see that, while the runs with the highest Mach numbers 
tend to be those with the highest overall final velocities, 
there are not necessarily the runs with the highest value 
of '^cioud/'^hot 5 meaning that the decrease of R_\_ often has 
a bigger effect on the final velocity of the cloud than the 
increase in the overall cloud lifetime. 

To study this evolution in greater detail, we normal¬ 
ized the cloud velocity in each case by the expectations 
from eq. ( 12 ) in the absence of changes in Av and R±, 
computing Ucloud(^) /'^estimate (^) with 


'^estimate (^) — 




(18) 


These ratios are plotted as a function of Mach number in 
Figurej^at each of the characteristic times defined above: 
^90 5 ^75 5 ^50 5 and ^ 25 , each plotted with a different sym¬ 
bol. This comparison shows that while Udoud/'^estimate 
varies significantly between runs, these differences, like 
the timescales in §4.1, are functions almost exclusively 
of Mhot- Furthermore, the dependence of Udoud/'^estimate 
is also only weakly dependent on the time at which the 
velocities are measured: with the earlier times being well 
described by Udoud/^estimate — 1.8(1 T Tfhot) and 
the late times being well described by Udoud/'^estimate = 
1.8(1 +Mhot)-'. 

These scalings can be understood in terms of changes 
in R±^ which means that a density increase by a fac¬ 
tor ^ 1 + Mhot will have a different effect on the cloud 
evolution depending on the direction of the contraction. 


If the contraction is purely in the direction of the flow, 
then R± will remain unchanged, and the contraction will 
not affect Udoud/'^estimate- On the other hand, a uniform 
contraction in density by a factor of 1 + M will result in 

^_l_ (X Tf ! such that Udoud/^estimate (^ T Tf) ^ . 

Finally, a contraction purely perpendicular to the flow 
will result in Udoud/'^estimate OC (1 + Mhot)-'. 

Both the Udoud/^estimate — 1.8(1 + Tfhot) and the 
'^cioud /'^estimate = 1 . 8(1 + M)“^ modcls plotted in Figure 
[^exceed 1 in the case in which Mhot < 1 - This means that 
they correspond to a increase in R^. In fact, a moderate 
lateral expansion, due to the initial shock that passes 
through the cloud, is visible in the low Mach number slice 
plots shown in Figure Comparing these slices with the 
higher Mhot runs in Figures and we can also see that 
the compression of the cloud at higher Mach number is 
somewhat biased to the direction to the flow at early 
times, resulting in a Udoud/^estimate scaling between (1 + 
Mhot)-^^^ and (l+Mhot)-'- At late times, the expansion 
of the cloud in the streamwise direction means that the 
compression is almost completely perpendicular to the 
flow direction, resulting in a scaling that goes as ~ (1 + 
Mhot)-'. 

In Figure we combine our fits for the mass loss and 
Udoud(^) to obtain 


Udoud(t) — 


3^hot[«(l + Mhot)'/"] 1.8 

(1 + Mhot)^ 


(19) 


where as in eq. a = 1.75, 2.5, 4, and 6 , at tgo, ^ 75 , tso, 
and ^ 25 , respectively and P = 0.8 or 1. This expression 
provides a good fit to all our data, with /3 = 0.8 fitting 
better at tgo and ^75 and P = I fitting better at late times. 
Together, eqs. ( p!4| ) and ( [T^ provide a good summary of 
our full set of results foiMhe evolution of clouds in the 
hydrodynamic case in which thermal conduction is small. 


4.3. Distances and Implications for Circumgalactic 
Obsevations 

Although conduction and magnetic fields will both 
play important roles in full models of cold cloud evo¬ 
lution, we can nevertheless make preliminary compar¬ 
isons between our results and current observations. From 


eq. (19), we can see that the cloud velocity goes as 


'^cioud ~ 1.4auhotXo ^ (1 + M)^-^“^, with a = 6 when 
the clouds are down to 25% of their original mass, rapidly 
shedding their remaining gas. The immediate implica¬ 
tion is that if xo ^ 100 no clouds will ever be accelerated 
to the hot wind speed before disruption. Furthermore, 
typical values of Udoud reach maximum values of O.Suhot 
or substantially less. 

Eq. (19), also shows that the cloud velocity at a given 


stage in its evolution depends almost completely on the 
velocity of the hot wind and the initial density contrast 
between the cloud and the exterior medium, with the 
largest velocities arising when Uhot is large and xo is 
small. To reach the ^ 200 km/s outflow velocities seen 
around large, rapidly star forming galaxies, {e.g. Heck¬ 
man et al 2000; Pettini 2001; Rupke et al 2002; Martin 
2005, 2012) thus requires that the clouds enter the hot 
medium outside of the driving region, where the wind 
velocities are the largest and the pre-shock ambient pres- 
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Fig. 6 .— Slices of pressure in units of ergs cm“^ (top row), density in units of g cm“^ (center row), and velocity in units of km s“^ 
(bottom row), taken at the midplane at six cloud crushing times. The slices are taken from runs Mlv860 (first column), M3.6v3000 
(second column), M6.2v3000 (third column), and M11.4v3000 (fourth column). All lengths are in kpc. 
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sures are the smallest. 

Integrating our fit to t;cioud out to the disruption times 
measured above, we can also predict the distance that 
the cloud will travel at different points in its evolution. 
This gives 


4ioud = 0.65a2i?^i„„d(i + 


( 20 ) 


with a and (3 defined as above. In Figure we com¬ 
pare this estimate with our simulation results. While the 
simulated distances show a bit more scatter than our ve¬ 


locity and timescale measurements, eq. (20) nevertheless 


gives a reasonable fit to the data, with /3 = 0.8 providing 
a good description at the earliest times, and /3 = 0.9, 
which averages over the velocity evolution seen at early 
and late times, providing a better fit at late times. As 
both of these [3 values are very near one, the distance the 
cloud travels as a function of mass loss remans roughly 
constant over a large range of external conditions. 


This means that, at least in the hydrodynamic case 
with no conduction, the distance traveled by the clouds 
depend almost exclusively only the initial cloud radius, 
with the clouds becoming completely disrupted as they 
move beyond ~ 40i?cioud- These distances are large 
enough that, given initial radii « 100 parsecs, cold clouds 
can travel from near the edge of the driving region to the 
few kpc distances at which they are typically observed 
in absorption against the starbursting host galaxy {e.g. 
Heckman et al 2000; Pettini 2001; Soto & Martin 2012). 

On the other hand, it is much more difficult for clouds 
to travel ~ 100 kpc distances, as probed in nearby 
galaxies using absorption lines measured in background 
quasars and galaxies (Bergeron 1986; Lanzetta & Bowen 
1992; Steidel et al 1994, 2002, 2010; Zibetti et al 2007; 
Kacprzak et al 2008; Chen et al 2010; Tumlison et al 
2013; Werk et al 2013, 2014; Peeples et al 2014; Turner 
et al 2014). Our results show that, regardless of the 
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Fig. 7. — Velocity evolution of the cloud as a function of Mach number and velocity. In each panel the solid (blue) line shows the mass 
weighted average velocity of the mass at or above 1/3 the original density of the cloud As in figure]^ in all panels the time is given in 
units of the cloud crushing time, and the length of each panel is also proportional to t/tcc- 


structure of the wind, such cold clouds would have to 
be the sizes of entire galaxies to travel to such large dis¬ 
tances without being disrupted. Thus, unless conduction 
and magnetic effects are able to preserve clouds for much 
longer than seen our current simulations, the ~ lO^K 
gas observed at ~ 100 kpc around galaxies can not be 
directly associated with lO^K ejected material. 

4.4. Resolution Effects 

As a test of resolution effects, we carried out two ad¬ 
ditional Mhot = 3.5 and Uhot = 1700 km/s runs: a 
high-resolution run with a maximum resolution of Ax = 
7?cioud/128, and a low-resolution run with a maximum 
resolution of Ax = i?cioud/32, as compared to our fiducial 
Ax = i?cioud/b4 resolution. Note that because we force 
all our runs to maintain the maximum resolution in the 
full volume surrounding the cloud, the high-resolution 
run maintains twice the resolution of the fiducial runs in 


all important regions of the calculation, while the low- 
resolution run has half the resolution of the fiducial runs 
throughout the most important regions. Our choice to 
force the regions around the clouds to the maximum res¬ 
olution also means the three runs vary dramatically in 
the computer time needed to complete them. To reach 
t = 12tcc, the low-resolution run took only « 3k CPU 
hours, the fiducial run took 24k CPU hours, and the 
high-resolution run took ^ 260k CPU hours. 

The upper left panel of Figure pT] shows the evolution 
of Ti/ 3 , the fraction of total mass at or above 1/3 the 
original cloud density, in each of these runs. The agree¬ 
ment between the fiducial and the high-resolution runs 
is excellent, with Fi/^{t) being almost indistinguishable 
for most of the evolution. On the other hand, the low- 
resolution mass evolution is significantly different, with 
the cloud taking much longer to be disrupted. This high¬ 
lights the importance of maintaining adequate resolution 
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Mach Number 

Fig. 8 .— Velocity evolution of the cloud normalized by 

'^estimate (^) = f^e expected value if the velocity of the 

4Xo ^cc 

exterior medium relative to the cloud, the radius of the cloud, and 
the density contrast between the cloud and the exterior medium 
all remain fixed. The points are velocities at tgo (violet asterisks), 
^75 (cyan pluses), (green diamonds), ^25 (black triangles). The 
blue solid curve shows 1.8(1 + Mhot)~^ and the red dashed curve 
shows 1.8(1+ Mhot)“°-®- 
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Fig. 9.— Velocity evolution of the cloud at times tgo (top 
left), ^75 (top right), ^50 (bottom left), ^25 (bottom right), 
normalized by Uestimate(^) = ^^ 2 ^^ evaluated at the estimated 

values of these times as given by equation ( |14| >. As in Figure 
the violet points are the simulation results, the blue solid 
curve shows 1.8(1 + Mhot)~^ and the red dashed curve shows 
1.8(1 + Mhot)~°'^- The more shallow curve provides the best fit 
at the earliest times, but at late times the results are better fit 
with the stronger scaling with Mach number. 


when studying the interactions described here. 

The evolution of the cloud velocity, shown in the top 
center panel of Figure is also very similar between the 
fiducial and the high-resolution runs, although, in this 
case, there are minor differences. While r’cioud(^) matches 
almost exactly between the two runs during the initial 
stages, at later times the cloud in the higher resolution 
run reaches slightly lower velocities than the cloud in 
the fiducial run. These differences, in turn, appear to 
be related to the evolution of material stripped from the 
core of the cloud. Finally, the late time history of the low- 
resolution run is much different than both of these runs. 



Fig. 10.— Distance traveled by the clouds as a function of the 
disrupted mass. The distances scale with cloud radius and thus 
are expressed in units i^cloud? which was 100 pa rsecs in our current 
simulations. The data is compared to eq. ( |20| ) with (3 = 0.8 (red 
dashed curves) and /3 = 0.9 (blue dot-daslied curves), with 0 .= 
1.75, 2.5, 4, and 6. 


reaching larger final velocities due to the underestimate 
of cloud disruption that occur in this simulation. 

In order to quantify the evolution of stripped material 
in each of these runs, we calculated the ‘mixing-fraction’ 
described in Xu & Stone (1995) and Orlando et al (2005), 
equation (19). This is defined as 

/mix ■ — - [ dVC,ioudP, ( 21 ) 

mcloud,0 Jy(0.1<Ccioud<0.9) 

where mcioud,o is the original mass of the cloud, and the 
integral is computed over all zones in the simulation in 
which the mass fraction of the tracer field Cdoud is be¬ 
tween 0.1 and 0.9. In this case, the fiducial and the high- 
resolution runs track each other well initially, but di¬ 
verge strongly after t ~ btcc , with the high-resolution run 
showing a larger mixed fraction than the fiducial case. 
The differences between the fiducial and low-resolution 
run are even stronger, with the low-resolution /mix ~ 0 
before 5tcc and about half of the fiducial /mix for most 
of the subsequent evolution. 

The lower panels of Figure show rendered images 
of the three runs at 5.8tcc, the time at which 75% of the 
mass is stripped away in the fiducial and high-resolution 
cases. Here we see that the distribution of the stripped 
material is significantly different between all three runs. 
In the fiducial run, the stripped material appears as spray 
of dense knots that maintain significant density contrasts 
out to large distances. On the other hand, the high- 
esolution run displays a smaller, wispier distribution of 
stripped material, consistent with better resolved mixing 
between the two media. Finally the low-resolution run 
shows most of the stripped mass collected up in four 
clumps arranged in an X pattern than is aligned with 
the diagonals of the simulation grid. 

These variations in morphology help to explain the dif¬ 
ferences seen in the evolution plots. The under-resolution 
of mixing in the Ax = i?cioud/32 run means that the 
retention of gas by the cloud is overestimated, leading 
to high mass fractions and low /mix values throughout 
the simulation. Furthermore, because the cloud remains 
coherent for longer times, r’cioud(^) is the higher in this 
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Fig. 11.— Comparison between M^ot = 3.5, i;hot = 1700 km/s, runs with maximum resolutions of i^cloud/32, i^cloud/64, and i^cloud/128. 
Top Left: Mass evolution in the fiducial run (blue solid lines) versus the high-resolution run (purple dashed lines) and the low-resolution 
run (red dotted lines). Top Center: Velocity evolution in the three runs, with lines as the top left panel. Top Right: Evolution of the 
mixed fraction, /mix 5 in each of the three runs. Bottom: Rendered plots of log density from the runs at a time ^75 = h.Stcc, with two 
surfaces plotted per decade over the range used in the slice plots above: log^^Q p/(gcm“^) = —25 to log^^Q p/(gcm“^) = —22.5. The plots 
are constructed with the camera pointed at the center of the cloud, viewing from a position {x^y^z) = (—3.3, —0.2, —2.0)Rcloud5 nnd were 
generated using the yt toolkit (Turk et al. 2011; http://yt-project.org). 


simulation than in the higher resolution runs. In the 
Ax = i?cioud/64, on the other hand, the cloud evolution 
is tracked well, but the stripped gas mixing is underesti¬ 
mated. Thus the more clumpy distribution in the fiducial 
case includes some fast moving material with densities 
above 1/3 of the initial cloud density, and this mate¬ 
rial lasts longer than it does in the high-resolution case, 
slightly increasing the measured cloud velocity. Finally, 
the evolution of /mix in the high-resolution run shows 
a larger fraction of mixed material that the other runs, 
as more cells are made up of wispy, 0.1 < Cdoud <0.9 
material. 

Computing mixing fractions for the Ax = i?cioud/b4 
with other Mach numbers and streaming velocities, we 
find that /mix varies from as little as 1 % to as much 
as 30%. In general, the higher the Mach number and 
the lower streaming velocity, the higher /mix, but these 
trends are noisy and uncertain. Thus while the Ax = 
i?cioud/b4 resolution adopted in our fiducial runs appears 
to be sufficient to measure the mass loss and velocity evo¬ 
lution in which we are most interested here, even higher 
resolutions are required to accurately track the structure 
of stripped material, which mixes into the hot medium 
at significant distances from the central cloud {e.g. Mari- 
nacci et al 2010; 2011; Kwak et al 2011; Henley et al 
2012 ). 

5. CONCLUSIONS 

Galaxy outflows are known to play a key role in the 
history of galaxy formation, but studies of their proper¬ 
ties are limited by the fact that the cold-cloud material 
that is easiest to observe is also the most difficult to un¬ 
derstand. As a first step in overcoming this limitation. 


we have conducted a suite of adaptive-mesh refinement 
simulations that include radiative cooling and track the 
evolution of such clouds under the full range of condi¬ 
tions relevant for galaxy outflows. By adopting large 
simulation volumes, carefully changing frame, and main¬ 
taining the highest resolution in the regions in which it is 
most needed, we have been able to accurately track these 
clouds for much longer times than previous simulations. 

This is particularly important for highly supersonic 
flows, in which the disruption of clouds by the KH in¬ 
stability is strongly suppressed. We find that cloud dis¬ 
ruption occurs at much later times in this case than 
in subsonic interactions, but still it occurs much sooner 
than expected from laboratory experiments of supersonic 
shear layers. Rather than finding lifetimes oc tccM~^ as 
expected from eq. ( 11 ), we find that the times at which 
that cloud loses 90%of its mass (tgo), 75% of its mass 
(^ 75 ), 50% of its mass (^ 59 ), tso and 25% of its mass (^ 25 ) 
are given by 

t = atcc \/l + Afhot, ( 22 ) 


where a= 1.75, 2.5, 4, and 6 , at tgo ^75 ^50 and ^ 25 , 
respectively. 

The reason for this scaling lies in the structure is the 
structure of the Mach cone that forms around the cloud. 
While the normal shock at the head of the cloud increases 
the pressure there by a factor 1 + M^, the downstream 
pressure rises only by a factor 1 + Mhot is it experiences 
only an oblique shock. Thus the cloud becomes com¬ 
pressed in the direction normal to the flow, but the pres¬ 
sure gradient in the streamwise direction causes it to ex¬ 
pand in this direction by a factor oc + Mhot)]- 

This in turn leads to a “cometary” appearance at inter- 
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mediate times and fi nally disru pts the cloud completely 
on a timescale oc tcc\/l + ^hot- 
We also find that the cloud velocity is primarily set 
by the momentum imparted by the impinging flow, as 
described by eq. ([T^. Our simulations show that devi¬ 
ations from this estimate depend almost solely on the 
Mach number. We find a fit for the cloud velocity as a 
function of the exterior Mach number given by 


'^cloud(^) — 


St’hot 1.8 




(l + Mhot)^’ 


(23) 


where again a = 1.75, 2.5, 4, and 6 , at tgo, ^ 75 , ^ 50 , and 
^ 25 , respectively, and where p = 0.8 fits best at tgo and 
^75 and /3 = 1 fits best at later times. These scalings 
can be understood in terms of changes in the effective 
radius of the cloud perpendicular to the direction of the 
flow. The (3 = 0.8 fit at early times represents cloud 
compressions that are somewhat biased to the direction 
perpendicular to the streamwise direction, while the the 
/3 = 1.0 fits at late time represents compressions that 
are almost completely perpendicular to the streamwise 
direction in which the cloud is stretched. 

Integrating our fit to 'T’cioud out to the measured dis¬ 
ruption times, we find ddoud = •b5(a^i?cioud(l 
with a and (3 defined as above. Because /d ~ 1, the dis¬ 
tance the cloud travels as a function of remaining mass 
fraction is roughly constant over a large range of exter¬ 
nal conditions, with the clouds becoming completely dis¬ 
rupted as they move beyond ddoud ~ 40i7cioud. These 
distances are large enough that, given initial cloud radii 
^ 100 parsecs, cold clouds can travel from near the edge 
of the driving region to distances of a few kpc, to be ob¬ 
served as absorption features superimposed on the spec¬ 
trum of there host galaxies. On the other hand, un¬ 
less conduction and magnetic effects are able to preserve 
clouds for much longer than seen our current simulations, 
~ lO^K gas observed at ~ 100 kpc distances around 
galaxies can not be directly associated with lO^K ejected 
material. 

In order to test resolution effects, we carried out runs 
with Mhot = 3.5 and Vhot = 1700 km/s, and maxi¬ 
mum resolution levels of Ax = i?cioud/32 and Ax = 
^cioud/128, as opposed to the fiducial values of Ax = 
^cioud/b4. We find that the Ax = i?cioud/32 results vary 
strongly from those of the other two runs, indicating that 
Ax = 7^cioud/32 is not adequate to obtain reliable results. 
On the other hand, the cloud mass and velocity evolution 
do not change significantly between the Ax = i?cioud/b4, 
and Ax = i?cioud/128 runs, although the mass fraction of 
the mixing layers between the fluids increases and their 
structure becomes more complex at very high resolution. 

The mass and velocity scalings found above can be im¬ 
plemented in numerical simulations of galactic outflows 
on larger scales. As the evolution of cold clouds in such 
simulations cannot be adequately modeled directly, sub¬ 
grid models that use our scalings for mass loss and cloud 
velocity could represent a promising path forward to to 
model the thermal evolution of the multi-phase material 
in outflows. However, before this can be carried out with 
confidence, we must first understand how thermal con¬ 
duction and magnetic fields impact the cloud and change 
its evolution relative to the hydrodynamic case with cool¬ 
ing. These physical effects will be simulated in forthcom¬ 


ing work. 


We would like to thank Romeel Dave, Paul E. Di- 
motakis, Timothy Heckman, Crystal Martin, Eve Os- 
triker, and Todd Thompson for helpful discussions. ES 
gratefully acknowledge the Simons Foundation for fund¬ 
ing the workshop Galactic Winds: Beyond Phenomenol¬ 
ogy which helped to inspire this work. He also grate¬ 
fully acknowledges Joanne Cohn, Eliott Quatert, and the 
UC Berkeley Theoretical Astronomy Center, Uros Sel- 
jak and the Lawrence Berkeley National Lab Cosmology 
group, and the organizers of the Gravity's Loyal Oppo¬ 
sition: The Physics of Star Formation Feedback at the 
Kavli Institute for Theoretical Physics, for hosting him 
during the period when much of this work was carried 
out. We would also like to acknowledge the Texas Ad¬ 
vanced Computing Center (TACC) at The University 
of Texas at Austin, and the Extreme Science and En¬ 
gineering Discovery Environment (XSEDE) for provid¬ 
ing HPC resources via grant TG-AST140004 that have 
contributed to the research results reported within this 
paper. Some simulations were also run on the JUROPA 
supercomputer at the Juelich Centre for Supercomput¬ 
ing under project HHB03. The FLASH code was de¬ 
veloped in part by the DOE-supported Alliances Cen¬ 
ter for Astrophysical Thermonuclear Flashes (ASC) at 
the University of Chicago. ES was supported by NSF 
grants ASTI 1-03608 and PHYl 1-25915, and NASA the¬ 
ory grant NNX09AD106. 


APPENDIX 


SHOCK JUMP CONDITIONS 

The are four distinct regions that describe the inter¬ 
action when the shock first hits the front of the cloud: 
the undisturbed cloud (region 1 ), the unshocked exterior 
medium (region 2 ), the region behind the shock that is 
transmitted through the cloud (region 3), and the region 
behind the shock that is reflected back into the exterior 
medium (region 4). These are illustrated in Figure 
(see also Silk &Solinger 1974; Hester et al 1994). Within 
the cloud, in the frame in which the cloud is initially 
stationary, we have: 

fh _ vt - V3 _ Mt - _ (71 - l)M‘^ + 2 

ps vt Mt (71 + 1 )M^ 

and 

^ ^ 271-^1 _ 7i - 1 
Pi 7i + 1 7i + 1 ’ 

where Vt is the velocity of the transmitted shock, Mt = 
and M 3 = xs/csp, and where we use ' to denote 
the ratio of the velocity behind the transmitted shock 
with the sound speed in front of the shock. In this same 
frame, within the exterior medium we have 


(Al) 

(A 2 ) 


P2 

Pa 

and 


V4-Vs _ M^- Mg _ (72 - 1)(M2 - + 2 

V2-Vg~ M2-Mg ~ (72 + 1)(M2 - M,)2 ’ 

(A3) 


P 2 


2^2{M2-Mgf 


72-1 

72 + 1’ 


(A4) 


72 + 1 
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Fig. 12.— Definition of coordinates for the shock. 


where Ms = Vslcs^ 2 , ^2 = '^ 2 /cs, 2 , and M 4 = v^lcs^ 2 , 
and again the ' denotes the ratio of the post shock ve¬ 
locity (in region 4) with the pre-shock sound speed (in 
region 2 ). We also have that ps = P 4 , 'T’s = 1 ^ 4 , and, 
because the pre-shock exterior material the cloud are as¬ 
sumed to be initial pressure equilibrium, ps/m = ^ 4 /^ 2 - 
Th ese relations allow us to equate eq. ( |A2[ ) and eq. 
(A4) to relate the incoming flow to shock passing through 
the cloud for arbitrary choices of 71 and 72 . If 71 = 72 = 
7, this becomes particularly simple. In this case we have 


that M 2 — Ms = Mf and so from eqs. (Al) and (A3), 
Pi/P3 = P 2 / pAi such that psip 4 = Xo- This means 


M, - M' = M' - M, = M' - M 2 + M,, 


(A5) 


and M 2 = M 3(1 -h Xo We then rearrange eq. (Al) 
by multiplying both sides by (7 + 1 )M| and plugging in 


for M 3 to obtain 


M/-l = M2Mt, 
7+1 

■ 


(A 6 ) 


where M 2 = M 2 - '^~ 7 i /2 ■ Solving for Mj, this gives 

2(l+Xo ) 


Mt = M2 


l + {l+AMpf/‘^ 


(A7) 


such that Vt ~ 1 ^ 2 /Xq_ _ the high Mach number limit. 

Combined with eq. ( |A 2 [ ) this gives the post-shock pres¬ 
sure in the cloud in terms of the original pressure. Tech¬ 
nically, these relations apply only when the Mach number 
of the incoming flow is 1 or above, but in practice extend¬ 
ing this to lower velocities gives reasonable estimates. 
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